# -*- coding: utf-8 -*-
"""
Created on Fri Mar 15 15:28:27 2024

@author: 86136
"""

import pandas as pd
import numpy as np
from matplotlib import pyplot as plt
from matplotlib.pylab import mpl
import tensorflow as tf
from sklearn.preprocessing import MinMaxScaler
from keras.layers import LeakyReLU
from keras.callbacks import ReduceLROnPlateau
from tcn import TCN, tcn_full_summary
from keras.callbacks import LearningRateScheduler
from keras.callbacks import EarlyStopping
from tensorflow.keras import Input, Model,Sequential

mpl.rcParams['font.sans-serif'] = ['SimHei']   #显示中文
mpl.rcParams['axes.unicode_minus']=False       #显示负号

#取数据
# 取数据
# 读取CSV文件，并将日期列解析为日期时间格式
data = pd.read_csv('文件路径')
data['Time'] = pd.to_datetime(data['Time'], format='IS08601')

# 提取年份和月份，然后重新构建日期格式列
data['date'] = data['Time'].dt.strftime('%Y/%m')

# 设置日期时间格式列为索引
datanew = data.set_index('date')

# 创建 Series，使用日期格式列作为索引
#series = pd.Series(datanew[['NJ', 'HF', 'SH']].values.flatten(), index=datanew.index)
series = pd.Series(datanew['QX'].values, index=datanew.index)

print(series)

#滞后扩充数据
dataframe1 = pd.DataFrame()
num_hour =24
for i in range(num_hour,0,-1):
    dataframe1['t-'+str(i)] = series.shift(i)
dataframe1['t'] = series.values
dataframe3=dataframe1.dropna()
dataframe3.index=range(len(dataframe3))

#二折划分数据并标准化
# pot=int(len(dataframe3)*0.8)
pd.DataFrame(np.random.shuffle(dataframe3.values))  #shuffle
pot=len(dataframe3)-40
train=dataframe3[:pot]
test=dataframe3[pot:]
scaler = MinMaxScaler(feature_range=(0, 1)).fit(train)
#s caler = preprocessing.StandardScaler().fit(train)
train_norm=pd.DataFrame(scaler.fit_transform(train))
test_norm=pd.DataFrame(scaler.transform(test))
test_norm.shape,train_norm.shape

X_train=train_norm.iloc[:, :-1]
X_test=test_norm.iloc[:, :-1]
Y_train=train_norm.iloc[:, -1:]
Y_test=test_norm.iloc[:, -1:]

#转换为3维数据 [samples, timesteps, features]
source_x_train = X_train
source_x_test = X_test
X_train = X_train.values.reshape([X_train.shape[0],1,X_train.shape[1]]) #从(909,16)-->(909,1,16)
X_test = X_test.values.reshape([X_test.shape[0],1,X_test.shape[1]])  #从(12,16)-->(12,1,16)
Y_train = Y_train.values
Y_test = Y_test.values

#动态调整学习率与提前终止函数
# 动态调整学习率与提前终止函数
#def scheduler(epoch, lr):
    # 每隔50个epoch，学习率减小为原来的1/10
    #if epoch % 50 == 0 and epoch != 0:
        #if lr > 1e-5:
           # lr *= 0.1
            #print("lr changed to {}".format(lr))
   # return lr

#reduce_lr = LearningRateScheduler(scheduler)
early_stopping = EarlyStopping(monitor='loss',
                               patience=20,
                               min_delta=1e-5,
                               mode='auto',
                               restore_best_weights=False,
                               verbose=2)

# 构造TCN模型
batch_size = None
timesteps = X_train.shape[1]
input_dim = X_train.shape[2] # 输入维数
tcn = Sequential()
input_layer = Input(batch_shape=(batch_size,timesteps,input_dim))
tcn.add(input_layer)
tcn.add(TCN(nb_filters=64,  # 在卷积层中使用的过滤器数。可以是列表。
        kernel_size=3,  # 在每个卷积层中使用的内核大小。
        nb_stacks=2,   # 要使用的残差块的堆栈数。
        dilations=[2 ** i for i in range(6)], # 扩张列表。示例为：[1、2、4、8、16、32、64]。
        # 用于卷积层中的填充,值为'causal' 或'same'。
        # “causal”将产生因果（膨胀的）卷积，即output[t]不依赖于input[t+1：]。当对不能违反时间顺序的时序信号建模时有用。
        # “same”代表保留边界处的卷积结果，通常会导致输出shape与输入shape相同。
        padding='causal',
        use_skip_connections=True, # 是否要添加从输入到每个残差块的跳过连接。
        dropout_rate=0.2, # 在0到1之间浮动。要下降的输入单位的分数。
        return_sequences=False,# 是返回输出序列中的最后一个输出还是完整序列。
        activation = LeakyReLU(negative_slope=0.01),
        #activation='relu', # 残差块中使用的激活函数 o = Activation(x + F(x)).
        kernel_initializer='he_normal', # 内核权重矩阵（Conv1D）的初始化程序。
        use_batch_norm=True, # 是否在残差层中使用批处理规范化。
        use_layer_norm=True, # 是否在残差层中使用层归一化。
        name='tcn' # 使用多个TCN时，要使用唯一的名称
        ))
tcn.add(tf.keras.layers.Dense(64))
tcn.add(tf.keras.layers.LeakyReLU(negative_slope=0.1))
tcn.add(tf.keras.layers.Dense(32))
tcn.add(tf.keras.layers.LeakyReLU(negative_slope=0.1))
tcn.add(tf.keras.layers.Dense(16))
tcn.add(tf.keras.layers.LeakyReLU(negative_slope=0.1))
tcn.add(tf.keras.layers.Dense(1))
tcn.add(tf.keras.layers.LeakyReLU(negative_slope=0.1))

reduce_lr = ReduceLROnPlateau(monitor='val_loss', factor=0.2,
                              patience=5, min_lr=0.001)

tcn.compile('adam', loss='mse', metrics=['accuracy'])
tcn.summary()

#训练
# 在model.fit()中添加callbacks=[reduce_lr]
history = tcn.fit(X_train, Y_train, epochs=5, batch_size=20, callbacks=[reduce_lr])
history.history.keys() #查看history中存储了哪些参数
plt.plot(history.epoch,history.history.get('loss')) # 画出随着epoch增大loss的变化图
# plt.plot(history.epoch,history.history.get('acc')) # 画出随着epoch增大准确率的变化图

#预测
predict = tcn.predict(X_test)
real_predict = scaler.inverse_transform(np.concatenate((source_x_test,predict),axis=1))
real_y = scaler.inverse_transform(np.concatenate((source_x_test,Y_test),axis=1))
real_predict = real_predict[:, -1]
real_y = real_y[:, -1]
pearson_corr = np.corrcoef(real_y, real_predict)[0, 1]
print(real_predict)
print("皮尔逊相关系数：", pearson_corr)
#误差评估
plt.figure(figsize=(15,6))
bwith = 0.75 # 边框宽度设置为2
ax = plt.gca()# 获取边框
ax.spines['bottom'].set_linewidth(bwith)
ax.spines['left'].set_linewidth(bwith)
ax.spines['top'].set_linewidth(bwith)
ax.spines['right'].set_linewidth(bwith)
plt.plot(real_predict, label='real_predict')
plt.plot(real_y, label='real_y')
plt.plot(real_y*(1+0.15), label='15%上限', linestyle='--', color='green')
# plt.plot(real_y*(1+0.1),label='10%上限',linestyle='--')
# plt.plot(real_y*(1-0.1),label='10%下限',linestyle='--')
plt.plot(real_y*(1-0.15), label='15%下限', linestyle='--', color='green')
plt.fill_between(range(0, 40), real_y*(1+0.15), real_y*(1-0.15), color='gray',alpha=0.2)
plt.legend()
plt.show()
